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Abstract. We revisit in this paper the location of the various components observed in the AGN of NGC 1068. Discrepancies 
between previously published studies are explained, and a new measurement for the absolute location of the K-band emission 
peak is provided. I t is found to be consi stent with the position of the central engine as derived bv iGallimore et all ll997l) . 
Capett i et alJll997l) and Kishimoto N 19991) . A series of map overlays is then presented and discussed. Model predictions of dusty 
tori show that the nuclear unresolved NIR-MIR emission is compatible with a broad range of models: the nuclear SED alone 
does not strongly constrain the torus geometry, while placing reasonable constraints on its size and thickness. The extended 
MIR emission observed within the ionizing cone is shown to be well explained by the presence of optically thick dust clouds 
exposed to the central engine radiation and having a small covering factor. Conversely, a distribution of diffuse dust particles 
within the ionizing cone is discarded. A simple model for the H2 and CO emission observed perpendicularly to the axis of 
the ionizing cone is proposed. We show that a slight tilt between the molecular disc and the Compton thick central absorber 
naturally reproduces the observed distribution of Hi of CO emission. 
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1. Introduction 

Efforts in understanding the physics at work in active galactic 
nuclei (AGN) started in the early seventie s, long after the first 
AGN or quasar had been discovered dSevfertl Il943t Isch midt . 
Il959l) . Great progress has been made over the last decade, 
thanks to improved observing tools and modeling explorations. 
On the side of observations, one can acknowledge the impact 
of better spatial and spectral resolutions in the near-infrared 
(NIR), mid-infrared (MIR) and radio wavelengths, as well as 
the access to high energy photons (UV, X-rays). On the side of 
modeling, the availability of complex codes taking into account 
the large number of physical processes taking place in AGN, 
and their inter-relations, has permitted to bridge in a more re- 
alistic fashion obser vations and model predictions (see for a 
review Eroiiklll999l) . 

It is now commonly accepted that the main source of en- 
ergy in AGN is of gravitational origin. The various components 
building up an AGN and its surroundings have also been iden- 
tified. However, their geometrical arrangement and their links 
remain to be understood in more detail, with the hope that this 
will bring also some clues about the formation process of such 
systems. The commonly accepted AGN model involves a cen- 
tral engine (black hole plus accretion disc, hereafter referred 
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to as CE) surrounded by sets of dense matter clouds (the so- 
called broad-line and narrow-line regions, hereafter BLR and 
NLR) and by a dusty/molecular disc-like and thick structure 
(hereafter the torus) which funnels the emission of high energy 
photons and particles along privileged directions (the ionizing 
cone, the radio jet axis). This in turn gives rise to viewing-angle 
dependent effects which are often invoked to explain the vari- 
ous brands of AGN discovered so far. And last, the AGN, in- 
cluding its close environment, are embedded in the bulge stellar 
population of the AGN host galaxy. 

The intrinsic physical scales of the different components 
mentioned above are quite different, with orders of magnitude 
such as: (a) CE (black hole plus accretion disc), 10 14 cm, (b) 
BLR, 10 I6 cm, (c) NLR, 10 20 cm, (d) torus, 10 18 to 10 20 cm, 
while the stellar system core hosting the AGN has a typical di- 
mension of 10 2l) cm (but can be more extended). On the other 
hand, the radio jet and structures can show a huge extension, 
in some cases on kpc scales. The different AGN components 
and the components contributing along the line-of-sight close 
to the AGN (such as the stellar population) can be disentangled 
in several ways: through their relative flux contributions which 
vary largely among AGN, possibly reflecting the AGN evolu- 
tionary state, and through their different angular sizes, accord- 
ing to the intrinsic sizes mentioned above. However, even in the 
close-by AGN NGC 1068, which lies at a distance of 14.4 Mpc, 
a torus component with intrinsic size 10 19 cm corresponds to an 
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apparent angular size of 0.04". So far, this is out of reach from 
an observational point of view, except in the radio range with 
interferometers, or using speckle techniques. 

Therefore, the comparison between observations and model 
predictions, necessary to probe model parameters, is limited by 
the angular resolution of available data sets. Conversely, one 
can focus on physical signatures specific to each component, 
some of which are listed hereafter. High energy photons and 
particles arise essentially from the CE (physical processes re- 
lated to matter accretion and emission at the surface of the ac- 
cretion disc), while the radio maser emission is thought to be 
related to the inner walls of the dusty/molecular torus. The flux 
variations of the CE provide information about the size and 
structure of the accretion disc. The impact of flux variations on 
the ionization state of the surrounding BLR clouds makes the 
BLR line emission respond with a time-delay which reflects the 
BLR size. The molecular line emission in the NIR to millimeter 
range, as well as the MIR continuum emission come predomi- 
nantly from the dusty/molecular component (torus or/and cold 
material in BLR/NLR clouds shielded from the intense UV 
radiation field). The contribution from the stellar population 
starts to fade off at wavelengths above 3.5 //m. In conclusion, a 
whole set of constraints can be used and this is why the multi- 
wavelength approach is a very powerful tool in AGN studies. 
Still, to assemble these pieces of information, we have to face 
the fact that they have been obtained at different angular reso- 
lutions and that they are most often lacking precise astrometric 
information. These are major sources of limitation which have 
indeed hampered so far a quantitative probe of model predic- 
tions. 

The preferred target on which one could solve the issue is 
of course the close-by AGN, NGC 1068. Because it is a type 2 
AGN, it contains all the components mentioned above; because 
of its proximity it offers two advantages, its brightness and one 
of the best intrinsic spatial resolution one can achieve at every 
wavelength. Thanks to these features, a very large number of 
observational studies have been performed on NGC 1068 since 
the early seventies (ADS returns 200 entries for refereed arti- 
cles containing NGC 1068 in their title, from 1970 to 2002). It 
is therefore tempting to try assembling major pieces of infor- 
mation obtained so far on this object, with the goal of locating 
its CE and analyzing the AGN surroundings. 

The first task is to register, on the basis of their absolute 
coordinates, high resolution maps available at different wave- 
lengths and locate the CE with respect to other AGN compo- 
nents. This is performed in Sect. |2 where we present a new 
determination of the absolute positioning of the 2.2 fj.m con- 
tinuum emission peak using ISAAC/VLT data in the H2 line 
emission and PdB/IRAM interferometric measurements in the 
CO(2-l) line and where we compile and discuss the few data 
sets with astrometric positioning available so far. After hav- 
ing secured the best possible registration, we compare maps 
obtained at different wavelengths and on different scales, in 
order to extract quantitative information to be compared with 
model predictions. In Sect. [5] we discuss the spectral energy 
distribution (SED) of the unresolved AGN component (intrin- 
sic size less than 8pc FWHM), while the origin of the MIR 
dust emission observed within the ionizing cone is investigated 



in Sect. @] Finally, we discuss and model the excitation of the 
molecular material located in the plane perpendicular to the 
ionizing cone axis (the torus plane) in Sect. [5] Concluding re- 
marks and prospective ideas are presented in Sect. [6] 

2. Relative positioning of the emission peaks at 
different wavelengths, location of the CE 

The precise location of the CE in the AGN of NGC 1068 is 
still subject to a hot debate. We review hereafter a number of 
published astrometric measurements at different wavelengths. 
Then, we present a new absolute positioning of the unresolved 
K-band continuum peak and we discuss the location of emis- 
sion peaks at different wavelengths with respect to the CE. 

2.1. The central engine 

It is commonly admitted that radio observations can directly 
probe the CE. The 6 cm MERLIN map of the AGN shows 
a linear radio structure made of four compac t components 
( Muxl ow et al ., 1996; lGallimore et al., 1996c a). One of these 
components (source SI) is now identified with little doubt 
as the CE, becaus e of its thermal inverte d spectrum between 
1.3 cm and 6 cm JGallimore et allfl996al) . Moreover, a disc of 
H2O and OH maser emission, thought to be associated with 
the inner walls of a molecular torus aroun d the CE, is ob- 
served with the VLBA at the same location jGallimore et all 
1996b). The ab solute astromet r ic prec ision of the MERLIN im- 
ages is 20 mas. Muxlow et al 1 11 19961) report the position of S 1 
at a = 02 h 42 m 40.7098 s , 6 = -00°00'47.938" (J2000). In the 
following, we consider that S 1 is the CE. 

2.2. Astrometric positionings 

The objective way to register any map to the 6 cm MERLIN 
image is by performing absolute astrometry. Any other method 
is subjective in the sense that it will require some physical as- 
sumption. 

Even though 0. 1" resolution HST UV-continuum and opti- 
cal NLR images have been pu blished before 1997 fevans et all 
ll99lHMacchetto et allfl994l) . they could not be precisely reg- 
istered since the absolute positioning of the HST images was 
only possible with a precis i on of +0.5". Based on a m ethod 
developed bv lLattanzi et all Jl997l) . ICapetti etall i 1997b report 
for the first time on an absolute registration of HST images with 
an ±0.08" precision. They find that the optical (F547M filter) 
continuum peak (hereafter referred to as HST F547M contin- 
uum) is located at a = 02 h 42 m 40.71 I s , 5 = -00°00'47.81" 
(J2000), that is 0.13"±0.1" North and 0.02"+0.1" East of SI. 
This corresponds to the peak of the NLR-cloud B. They also 
show that the UV (F218W) and F547M peaks are coincident 
within 25 mas. 

On the photogra phic 103a-O platest hey use to perform the 
absolute astrometry, Cape tti et al I II 19971) locate the visible pho- 
tocenter at a = 02 h 42 m 40.727 s , 6 = -00°00'47.69" (J2000). 
The difference between the locations of the HST F547M con- 
tinuum peak and the photographic plate visible photocenter is 
0.27". Such a large difference can arise from two effects. First 
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Fig. 1. Illustration of the seeing effect on the [OIII] image of the ionizing cone in NGC 1068. From left to right: HST [OIII] 
image, same image convolved with a 0.5" Gaussian, same image convolved with a 1.0" Gaussian. The "+" marks the position of 
the HST [OIII] peak and the "x" marks the position of the [OIII] peak in the images degraded by seeing effects. 



the F547M filter cuts out the prominent [OIII] /15007 emis- 
sion lines, while the ground-based 103a-O photographic im- 
age does sense these lines, the photocenter of which may be 
slightly displaced from the continuum photocenter. However, 
such a displacement of the [OIII] peak with respect to the 
F547M peak only occurs at the scale of individual NLR clouds. 
Indeed Macchetto et al. ill 994 show that the UV continuum 
emission peak in NLR-cloud B is shifted by 0.05" towards the 
North with respect to the [OIII] peak. Actually, the main rea- 
son for the observed shift of 0.27" is that the photographic im- 
age, which has a low spatial resolution (larger than 1"), mixes 
the different continuum and [OIII] emitting clouds. Indeed, one 
observes that the photographic plate photocenter is shifted to- 
wards the North-East with respect to the F547M peak. This 
is consistent with the overall shape of the NLR, which is 
extended towards the North-East. Considering the complex- 
ity of the NLR visible continuum and [OIII] emission, im- 
ages at different spatial resolutions ought to produce emission 
peaks at different locations. This effect must be considered 
carefully when interpreting the relative astrometric measure- 
ments between the emission peaks in the optical and in the 
IR, since the bul k of the NIR em ission is unresolved down 
to 0.03" (Wittkowskietalil 1998), while the optical emission 
is made of HST-resolved clouds over a region of ~0.5" xl". 
Fig.n illustrates t his effect by showing th e 0.1" HST [OIII] 
image taken from iMacchetto et alJ l l 19941) . and the same im- 
age after convolution with 0.5" and 1.0" FWHM Gaussian 
profiles, simulating seeing effects. The shifts of the position 
of the [OIII] peak on the images respectively degraded by an 
0.5" and 1.0" seeing effect are 0.35" and 0.45", with respect to 
the HST original image. This is of the same order of magnitude 
as the measured offset between the HST F547M peak and the 
103a-O visible photocenter, sensing the continuum (and [OIII] 
lines for the 103a- O pl ates). As a consequ ence, the results of 
iBraatz et all i 19931) and lMarco et all i 19971) which position the 



N-band peak and K-band peak at 0.3" South, 0.1" West from 
respectively the optical peak (from an optical TV camera) and 
the I-band peak should be considered with caution: they both 
positioned the MIR or the NIR peaks with respect to the visible 
peak of a low resolution image. In order to gain in precision, 
we discuss in the following the positioning of the NIR peak 
with respect to the optical continuum peak of NLR-cloud B 
only. 

2.3. Assumption-depending positionings 

Positioning of the emission peaks at different wavelengths with 
respect to the CE have also been performed under various phys- 
ical assumptions. 

2.3.1 . Positioning of the K-band peak 

The absolute astrometry of the K-band emission peak can be 
performed directly in two ways. 

The first method is to compare the high resolution CO emis- 
sion line map, for which absolute astrometry is available, and 
the H2 2.12 jjm emission line map, which appears to be very 
similar to the CO m ap and for which ther e is a relative pos ition 
of the K-band peak. lTacconi et alJ Jl997h and lSchinnerer et all 
(2000) observed at ~0.5" resolution the CO(2-l) line emission 
in the millimeter wavelength range with the IRAM interferom- 
eter. They produced a CO map of the central 3" region, find- 
ing that the CO emission is asymmetrically distributed, with 
a strong CO emission knot to the East of the CE. The CO 
emission peak corresponding to this knot is located at a — 
2 h 42 m 40.767 s and 5 = - 0°00'47.95" (J2000), wit h a preci- 
sion of +0.05". In parallel. [Galliano & Alloinl (l2"002) produced 
a map at similar resolution of the H2 2.12 fim emission derived 
from VLT/ISAAC observations. The CO and H2 maps corre- 
late very well (Fig.^e). Measurements on the data presented 
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in iGalliano & Alloinl J2002) show that the eastern H2 knot is 
located 0.95"±0.15" East and 0.1" ±0.15" North of the K- 
band continuum peak. The similarity of the two maps and 
the fact that, in the AGN context, CO and H? emission are 
expected to be coincident jMalonev et all Il996l) . allow us to 
assume that the CO and H2 peaks are at the same location 
within 0.1". The absolute coordinates of the K-band peak 
can therefore be derived: a = 2 h 42 m 40.702 s ± 0.012 s and 
6 = -00°00'48.04" ± 0.18" (J2000). This position of the K- 
band emission peak is consistent with the position of SI, while 
it is only marginally consistent with the positions of the HST 
F547M continuum and the [OIII] peaks. 

Another positionin g of the K-band peak has been per- 
formed by Thomp son et alJ ( 1200 ll) . Using the similarity of the 
NLR features seen in the HST/WFPC2 H a image and in the 
HST/NICMOS P a image and using the fact that both U a and 
P„ are Hydrogen recombination lines, and consequently should 
trace the same features, they have registered the position of the 
NIR continuum image with respect to the HST F547M con- 
tinuum image. They find the location of the K-band peak to 
be consistent with both the center of th e UV polarization pat- 
tern as determined by Kishi motol ill 9991 see Sect. 12. 3> . and the 
[OIII] peak (NLR-cloud B). This coincidence is claimed to be 
accurate to ±0.04". 

2.3.2. Positioning of the center of the polarization 
patterns in the UV and in the IR 

Capet ti et alJ<1995albh have presented HST imaging polarime- 
try of the NLR in the UV (F218W) and optical (F555W). 
The polarization vectors are found to be distributed in a cir- 
cular pattern, as expected in the case of scattered radiation 
from a point-like source, assumed to be the CE . After a cor- 
rection of their first result JCapetti et al they claimed 
to have determined the position of the center of polarization, 
at the location 0.45" ±0.015" South, 0.06±0.01 5" West of 
the op tical (F555W) continuum peak. However, iKishimotol 
(1999) noticed th at in the map of polarization position an- 
gles (PA) o f lCapettietalllll995albl) . clear deviations from the 
centrosymmetric pattern can be seen , although not discussed 
in that paper. Therefore, Kishimoto ( 1999) has extended this 
work, performing a very complete error budget, and has re- 
determined the location of the symmetry center of the po- 
larization pattern in the UV within a convincing error circle. 
He finds a locatio n for the CE wh ich is significantly different 
from that of Cape tti et a £l995b). The revised location is at 
a = 02 h 42 m 40.713 s ± 0.0006 s and <5 

(J2000). 

iLumsden et alJ Jl999h have performed NIR and MIR imag- 
ing polarimetry at good spatial resolution (~0.5" in the NIR). 
They detect NIR polarized emission in extended regions along 
PA~30° both towards the North and the South. In the optical 
the light from the southern part is absorbed, probably by matter 
in the plane of the host galaxy. They suspect that an additional 
mechanism, other than scattering, is contributing to the polar- 
ization at all wavebands close to the AGN. Therefore, they have 
masked the central 2" region and derived the scattering center 



-00°00'47.93" +014 „ 
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from the outer regions only. They find offsets between the po- 
larization pattern center and the flux centroid in the three NIR 
bands (J,H,K) ofrespectively0.il" ±0.18" , 0.09" ±0.21" and 
0.34" ±0.39" , offsets which are formally consistent with no 
offset. Unfortunately, they do not provide the PA of the offsets. 
If the PAs for the three offsets are different, then their measure- 
ments will actually be consistent with no offset, while on the 
contrary, if the three offsets were for example at PA~30°, then 
their resu lt would be consis tent with the UV polarization center 
found bv lKishimotol Jl999h 

Recently, liimpson et all J2002I) have published NIR (1.9- 
2.1 fim) polarimetry from HST/NICMOS data set. They claim 
that the center of the NIR polarization pattern does correspond 
to the NIR peak within ±0.08". 
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Fig. 2. Positions of the different peak, as described in Sect.|2] 
on the J2000 absolute coordinate reference frame. The errors 
represented by the boxes are given with respect to the absolute 
coordinate reference frame. The hatched rectangle marks the 
location of SI (CE). The dotted line rectangle corresponds to 
the location of the F547M continuum peak. The cross marks 
its center, while the diamond gives the relative position of the 
[OIII] peak (NLR-cloud B) with respect to the F547M peak; 
the size of the diamond gives the error of this relative measure- 
ment. The plain line rectangle corres ponds to the center of the 
UV polarization pattern center (Kishimoto, 1999). The dashed 
line rectangle represents our new measurement of the K-band 
peak location. 



2.4. Summary and conclusion 

We have compiled measured positions of the emission peaks 
at different wavelengths, of the centers of the polarization pat- 
terns in the UV and in the IR and of the positions of the radio 
sources in the region surrounding the AGN in NGC 1068. We 
find that the measurements converge to a unique picture, and 
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RA (J2000) 




dec (J2000) 


reference 


F547M peak 


02 h 42 m 40.711 s ± 


0.005 s 


-00°00'47.81" ±0.08" 


CaDetti et al. 0997) 


SI 


02 h 42 m 40.7098 s : 


± 0.001 s 


-00°00' 47.938" ± 0.02" 


Muxlow et al. (1996) 


K peak 


02 h 42 m 40.702 s ± 


0.012 s 


-00°00'48.04" ±0.18" 


this paper 


UV polarization 
pattern center 


02 h 42 m 40.7 13 s ± 


0.006 s 


-00°00' 47.93" +0M ", 

-0.11 


Kishimoto H999) 



Fig. 4. (flgure4.jpg) Map superpositions. Axis unit are arcsec, North is up and East to the left. Coordinates (0";0") give the 
location of the CE (see Sect.[2). F rom top to bottom, left to right: (a) 12"xl2" region around the CE of NGC 1068. T he gre y 
contours show the 6cm emission JGallimore et al.[[l996ch : the red dotted contours show the 20 fim image JAlloin et all [20001. 
(b) Zoom to the c entral 4"x4" region. Dotted red contours are the same contours as i n map (a); the so l id red contours show 
the 12.3jt/m map (jTomono et al., 2001) . (c) In green is displayed the [OIII] image by Ma cchetto et alJ (|l994); in black con- 
tours the 6 cm map b y Gallimore et al. ( 1996c); red contours are the same as in (b). (d) The dashed contours correspond to the 
Bry image by|Galliano & Alloin (2002 ); the green image and the black thin contours are as in (c). (e) The orange map shows 
the H? 2.12 ^m emission from lGalliano & AlloiiJll2002l) : the large scale contours show the 12 CO(2-l) map from S c himiererej^lj 
(2000); the small scale contours are the same as in map (c); the blue map displays the M-band image bv lMarco & Alloinl t2000) 
(f) The orange image is as in (e); the red contours as in (b); the dashed contours as in (d) 



Fig. 3. (figure3.jpg) Superposition of the [OIII] image of the 
NLR by Macchett o et all Jl994l) (grey scal e), MERLIN 6cm 
image (thick conto urs Muxl ow et all 1 19961) . and 12.7 yum im- 
age by B ock et all J2000I) in thin contours, following the posi- 
tions given in tabled The la bels are f rom Evans et al. (1991) 
for the NLR clouds and from lGallimore et all Jl996d) for the 
radio sources. Coordinates (0",0") are assigned to SI. 



tered by the dust along the axis of the torus, where the optical 
depth to the observer is smaller, and consequently, the center 
of light in the NIR may be displaced towards the scattering re- 
gion, away from the CE. As the wavelength increases in the 
MIR, the optical depth to the CE decreases, as does the scatter- 
ing efficiency, and the MIR light peak gets closer to the location 
oftheCE. 



that all the components are relatively positioned with a good 
precision. A sketch of the positions of the most important fea- 
tures is displayed in Fig.[2]while their absolute coordinates are 
listed in Table d 

iThompson et all ( l200lh found that the K-band peak is co- 
incident with the position of the UV polarization pattern center 
and the p eak of rOIHU5 07 in the NLR-cloud B, at +0.04" 
precision. Simps on" et all J2002I) found that, within +0.08", 
the NIR (1.9-2.1 //m) polarization pattern is centered on the 
NIR peak. These three measurements converge to position 
the NIR peak at the position of SI within ±0.1". We do 
not contemplate the offsets between visibl e peaks and N- and 
K-b and peaks m e asured respectively by Bra atz et all (11993) 
and iMarco et all Jl997t) . since these offsets are affected by 
the use of images at different spatial resolutions (Sect. 12. 2> . 
Regarding the UV polariza t ion pa ttern center, we consider only 
the analysis of iKishimotol Jl999h . which is the most precise. 
There is no absolute positioning of the MIR maps available. A 
very reasonable assumption, supported by dust y torus models 
JCranato fc n a ne S elll994lGT a nato et alltl997li. is to consider 
that the MIR emission peak appears closer to the CE than the 
NIR emission peak. The reason for this is that the optical depth 
in the NIR is larger than in the MIR and the scattering efficiency 
is much higher in the NIR than in the MIR. Therefore, the NIR 
photons that we see might well be photons which were scat- 



Fig-El displays a superposition of the radio 6cm contours 
and the MIR 12.7 /urn contours over an image of the NLR in the 
[OIII]/15007 line emission. The positions are set according to 
Sect. 12.41 The coincidence found between the NLR [OIII] emit- 
ting clouds and the MIR emission suggests that the same phys- 
ical clouds are responsible for both emissions. With this regis- 
tration, the radio source C is not coincident with the NLR-cloud 
B, as su ggested in the preferred registration of Galli more et all 
( 1996c). If one attempts to identify a particular NLR cloud re- 
sponsible for the deflection of the radio jet, it would more likely 
be NLR-cloud C. The fact that NLR-cloud B does not affect 
the trajectory of the jet suggests that it lies on the foreground 
(or background) of the jet. The anti-correlation between the 
[OIII]+MIR images on one hand and the radio map on the other 
hand, is striking. It looks as if the [OIII] and MIR light were 
emitted principally in regions surrounding a cavity carved by 
the radio jet and localized at the edges of the ionizing cone. The 
matter inside the cone may have been swept away by the jet, 
and pushed to the actual location of NLR-clou d G, 1 .7" f rom 
SI in th e direction of the deflected radio jet. Galli more et alJ 
Jl996d) originally proposed that the gaps in the NLR clouds 
may be caused by the jet, guessing an alignement between the 
MERLIN and the HST images, that turns out to be close to the 
actual alignement, as found in our analysis. 
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2.5. Map superpositions 

Considering the absolute and relative positionings discussed 
above, Fig. |4] presents a series of superpositions of maps at 
different wavelengths and on different scale. The caption de- 
scribes each map in detail. 

Map (a) shows the excellent correlation between the large 
scale 6cm map and the 20 fj.m image. At this resolution (0.6"), 
each 20 fim peak has a counterpart on the 6cm map. 

Map (b) reveals the good correlation existing between the 
low resolution (0.6") 20jum image and the high resolution 
(0.1") 12.3 jt/m map. This demonstrates the very good consis- 
tency of the two datasets. 

Map (c) shows the correlation between the [OIII] and the 
MIR emission and the anti-correlation between the [OIII] and 
the 6 cm radio emission on small scale. This is a strong evi- 
dence in favor of ionized clouds containing large amounts of 
cool dust. 

Map (d) indicates that, as expected, the Bry emission is 
correlated with the [OIII] emission. It also shows that the 
Bry emission is strongly enhanced at the location of NLR- 
cloud E. Its particular excitation, the fact that it is not a strong 
MIR emitter (see Fig.[5} and its location with respect to the 
radio jet suggest that NLR-cloud E is undergoing a different 
excitation than the other NLR-clouds. 

Map (e) highlights the very good correlation between the 
H2 2.12//m and the 12 CO(2-l) maps. Still one notices the in- 
verted brightness ratios in CO and H2 in the two western knots. 
On the contrary the eastern knot is very bright in both lines. 
This map also shows the spatial coincidence, at a 0.5" scale, 
between the radio emission and the M-band emission. One even 
notices that the M-band emission bends together with the radio 
emission, at about 0.3" North of the nucleus. 

Map (f) summarizes the locations of the H2 2.12/im and 
Bry line emission, and the 12.3 fim emission. 



3. The unresolved torus emission 

The observed IR SED of NGC 1068 and other AGN have 
been compared to predictions from radiative transfer models 
within dusty tori, in order to investigate the validity of the 




g structure (e.g. Pier & Krolik, 1993; Granato & Danese, 


Efstathiou & Rowan-Robinson, 


1995; 


Granato et al., 


Alonso-Herrero et alj. 


2001: iNenkovaetall 120021). A 



major problem is that it is still quite difficult to assess the 
truly nuclear SED to be compared with torus models, even 
for the best studied AGN, NGC 1068. In the NIR bands, even 
small aperture photometry is likely to include a non-negligible 
contribution from starlight. For this o bject, the best option so 
far is to combine the nuclear fluxes by Rouan et al. ( 1998J) in 
the J-,H-,K-bands (cor rected for star contamination) and by 
iMarco & Alloinl feOOOh in the L-,M-bands with the nuclear 
MIR fluxes provided in lTomono et alJ J200ll) since they are all 
derived from high spatial resolution images. The resulting nu- 
clear SED is shown on Fig. [5] and the flux values reported in 
Tabled 




Fig. 5. Examples of fits for tapered disc models. The model pre- 
diction is shown as a continuous line: from top to bottom, the 
model parameters corresponds to models 1,2 and 3 in Table[5] 
The observed SED is represented by diamonds (plus error bar) 
and corresponds to the fluxes displayed in Table|2] 



In the absence of precise physical ideas concerning the 
structure of the obscuring torus, several geometries - and as- 
sociated free-parameters - are plausible and indeed have been 
investigated in the papers quoted above (flared discs, tapered 
discs, cylinders, with and without substantial clumping etc.). 

An usually overlooked problem in this kind of study, is that 
dust optical properties around AGN are likely to be to some 
extent peculiar. As a matter of fact, Maio lino et alJ (2001a b) 
presented evidence for 'anomalous' properties of dust in AGN. 
Here the term 'anomalous' is relative to the standard prop- 
erties of dust grains thought to be responsible for the aver- 
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Table 2. IR flu xes for NGC 1068 The re ferences are the 
follow ing: R 98=dRouan et all 1 1 99 A M00=i lMarco & AlloinL 
2000), T01= (lTomono et all 1200 lb 
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age Mi lky Way extinct i on law an d cirrus emission. As pointed 
out by IMaiolino et al . (2001 albl) . it is not surprising that the 
properties of dust grains might be very different in the dense 
and extreme environment of an AGN. In particular, they in- 
voke a dust distribution biased in favor of large grains, rang- 
in g up to 1-10 um, mu ch larger than the usual cut at ~ Q.3um 
oflMathis etalJlll977h type models (e.g. lDraine & Leetfl984l 
Silv a et all [1998). This solution has been introduced as the 
best, among the many explored, to explain the well established 
decrease of E B -v/Nh and A V /N H ratios in the line-of-sight of 
AGN. It is worth noticing that the tendency for larger Ry = 
Av/(Ab - Ay) found in den se galactic reg i ons are commonly 
expla ined in a similar way dPrainei l200lt IMaiolino & Nattal 
2002, and references therein). 

Thus, we have performed SED fitting, adopting several dif- 
ferent geometries, but also investigated models in which the 
size distribution of grains extends to radii, a max , larger than the 
standard value. Of course, the uncertainty on optical proper- 
ties of dust further limits the already moderate success of this 
approach in constraining the geometry (see below). 

Due to the relatively long computing time required to pro- 
duce a model, the SED fitting is performed through compari- 
son with libraries of models. Each library consists of several 
hundreds of models belonging to a given "geometry class", 
in which typically 4-5 parameters are varied assigning to 
them 3-4 different values over a quite wide range. In par- 
ticular we have compared the NGC 1068 nuclear SED with 
both anisotropic flared discs and tapered d iscs (for definitions, 
see lEfstathio u & Rowan-Robinsoni I 19951) . The former geom- 
etry consists of a structure whose height above the equato- 
rial plane h increases linearly with the distance in the equato- 
I 1 1 ' 1 

rial plane r (Fig. la in Efstathiou & Rowan-Robinson, 1995). 

Therefore the dust-free region (ignoring here dust possibly 
surviving in the NLR clouds) is exactly conical, with half- 
opening angle ®h. We also introduce (at variance with respect 
to lEfstathiou & Rowan-Robinsoni 1 19951) a dependence of the 
dust density p, and therefore of the radial optical depth, on the 
polar angle O. The general form used for this dependence and 
that on the spherical coordinate r is: 



p(r, 0) oc r 13 exp(-a- cos 2 0) 

where a and /3 are model parameters. 

In tapered discs, h increases with r only up to a max value 
^max- In this case the dust density within the disc has been al- 
ways kept independent of ©. The two additional parameters 
used to fully characterize the models in both cases are the ratio 
between the outer and inner (i.e. sublimation) radii r max /r m i n , 
and the equatorial optical depth to the nucleus at 9.7yum, T e pj 

The parameter values used in the libraries are: 

- for both flared and tapered discs: 

O-max - 0.3,0.5, l.Oyum 

/? = -1,0,0.5,1,2 
T e ,9. 7 =0.1,0.3,1,3,10 

r max I K min 

= 30,100,300 

- for flared discs only: 

0/, = 20,30,40° 
0^ = 0,1,3,6 

- for tapered discs only: 

®h = 30° 
h max /r aax = 0.2,0.3,0.4 

The impacts of changing variou s parameters on t h e pre- 
dicte d SED have been e xplored in Grana ~& Danesel |[l994) 
and iGranato et alJ dl997l) . Here we recall only that, roughly 
speaking, r mflA /r m ,„ is related to the width of the IR bump, 
whilst T e ,9 .7 controls mainly the NIR slope of the SEDs, as ob- 
served from obscured directions, as well as its anisotropy. 

All possible combinations of parameter values are consid- 
ered in the libraries. For each corresponding model we com- 
pute the SED for viewing angles from 0° to 90° in steps of 10 
degrees. Then we search through these libraries which models 
provide a reasonable fit to the observed SED. The employed 
merit function is x 1 - We have assigned to each data point an 
error of 30%. This could well be an underestimate of the true 
uncertainty of the nuclear flux, especially in NIR bands, due 
to the difficulty of disentangling the starlight contribution. For 
instance, the nuclear K-band flux (corrected for stellar light) by 
Alonso-Herrero et al. (2001) in 3" aperture is greater than our 
value by as much as a factor 3. 

Indeed, the main result of our analysis is that several 
rather different combinations of geometry and parameters of 
the torus are perfectly compatible with the observed IR SED of 
NGC 1068 (examples of library fits are shown in Fig.|5]for the 
tapered models). In other words, as expected, the nuclear SED 
alone does not allow to put strong constraints on geometrical 
parameters. Several previous analysis have been too optimistic 
on this aspect, due essentially to an insufficient exploration of 
the parameter space. There is some tendency though, for the 
models with a max larger than the standard value to allow more 
combinations of the other parameters and still get an acceptable 
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fit. We also confirm that some constraint can be put on the torus 
extension and on the optical thickness, and we find that good 
fits are obtained both with tapered (fits shown in Fig.|5} and 
flared disc models, provide d the tori are m oderately thick, in 
the sense discussed by Gran ato et alJdl997h . The development 
of interferometric observations should, in the future, provide 
more observational constraints. 



Table 3. SED fitting parameters. The parameters meaning is 
given in the text, "t" refers to tapered discs, and "f" to flared 
discs. 
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0.5 


0.5 
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0.0 


0.5 


0.5 



4. The extended MIR emission associated with the 
ionizing cone. 

MIR high resolution images of NGC 1068 show elonga- 
tion perpendi cular l y to t h e putative d usty torus-like struc- 
ture JAlloin et all hoOfj bock et all l2000t iTomono et all 
l200ll) . This is interesting in two ways: (i) the ex- 
tension of the emission, of the order of 70 pc, ex- 
ceeds the estimated size of the torus, which is not ex- 
pected to be larger than 20 p c dGranato & DaneseL 1 19941 
lEfstathiou & Rowan-Robinsorl Il995t iGranato et all fl997l) . 
and (ii) it traces approximately the two ionization cones where 
dust grains are less likely to survive. Also, this emission is more 
prominent on the North-East side of the CE (which is believed 
to trace the cone closer to our line of sight) than on the oppo- 
site South- West side. The decline of the surface brightness with 
distance from the CE is relatively shallow and the observed "lo- 
cal" SEDs of the emission (see Fig.|6jl show interesting trends: 
in the North-East side, the SEDs along the axis have a steep 
decline towards shorter wavelengths, while in the South- West 
side, the SEDs appear to be flatter. 

4.1. Diffuse dust in the cone? 

Is it possible to explain the observations with diffuse dust dis- 
tributed in the cone? To investigate this hypothesis, we have run 
several sets of models with diffuse dust in the cone extending 
up to ~ 1", with different density and grain size distributions, 
and could not find a good match with the observations. Such 
models are affected by either, but usually both, of the follow- 
ing problems: (i) the emitted specific intensity decreases more 
steeply than the observed MIR maps show, and/or (ii) the local 
SEDs in the cones are too flat. The dust volume emissivity can 



be enhanced either by increasing the dust density or its tem- 
perature. The first possibility works only up to a certain point 
above which the dust in the cone becomes optically thick to 
optical-UV photons and, as a consequence, is less heated by 
the CE. We have tested with our model that this "turnover" is 
reached before the observed emission can be produced (see also 
the analytical estimate below). Moreover, an excess density of 
diffuse dust in the cone would obscure the CE and the BLR 
from any direction, and would also require a different mecha- 
nism than direct photo-ionization from the CE to produce the 
conical NLR. 

Let us perform as well an analytical calculation in order 
to evaluate when this "turnover" is reached. Let us assume a 
distribution of homogeneous diffuse dust in the cone, optically 
thin in IR (not necessarily in optical-UV). The observed spe- 
cific intensity I Y< \% at 18 fim is: 

7 v ,i8 = la V) i$B v ,n( T ) = T v,i8(0-B v ,i 8 (r) 

T y> i8(0 is the optical depth along the line of sight crossing the 
cone (1), ff v ,i8 the extinction coefficient and B v is the Planck 
function at the frequency v. Specific quantities are all evalu- 
ated at 18/im. Given the geometry, we can evaluate the dust 
temperature T assuming that the heating is only from the CE 
and that the distance from the central source, r, is almost con- 
stant along the line-of-sight. Also, r I. Thus we can identify 
the optical depth r v jg(/) with the optical depth to the central 
source T V; is(r). For graphite grains of typica l size, the equilib - 
rium dust temperature T can be estimated as (B arvainislll987l) : 



T = 1650 



(L, 



UV,4(> 



pc 



1/5.6 



exp -- 



Tuv(r) 
5.6 



Setting r — 70 pc (i.e. 1", where the observed 7 Vj i8 is 
Jy arcsec~ 2 ) and Luv = 0.02 x 10 46 erg s we have 
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T = 180 exp - 



Tuv(r) 
5.6 



On the other hand, for standard Galactic dust mixture, Tuv 
80 x t v is. This leads to: 



'v,18 



Tuv(r ) 
80 



fl„,i8 180 exp 



Tuv(r) 
5.6 



This function of Tuv has a maximum of ~ 10 Jy arcsec~ 2 (i.e. 
more or less the observed level) at tuv - 1, but as long as the 
dust is optically thin in the cones at optical-UV wavelengths - 
the condition for photoionization - the predicted 18 fim bright- 
ness is significantly lower than observed (for instance by a fac- 
tor 4 for Tuv - 0.1, or by a factor 100 for t V v — 0.01). 

To increase the dust temperature, it is necessary to invoke 
some additional heating source for the dust in the cone, such as 
anisotropic emission from the CE, young stars, or even s hock 
heating (see also discussion in I Villar-Martm et all l200ll and 
Sect. I4.3i . While all these are rather likely to occur, they are 
not satisfactory explanations in this respect: indeed the dust 
temperature in the cone, under the assumption of optically thin 
dust heated only by an isotropic central source, is already too 
high (~ 200-300 K at ~ 0.5", depending on the dimension and 
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Fig. 6. Background image: the 12.3 fim map from lTomono et all J200ll) . In each cell of the image, we show the observed local 
SED between 7 jjm and 20 //m. Each SED is the total SED of the corresponding cell. The units for the SEDs are Jy for the y-axis 
and yum for the x-axis, with the scales given in the insert on the top-right. 



composition of grains), to produce a SED declining as steeply 
as observed towards short wavelengths. A higher T, as could re- 
sult from an additional heating source, would worsen the prob- 
lem, independently of the nature of this source. 

Therefore, we are led to conclude that a distribution of dif- 
fuse dust within the ionizing cone cannot account for the MIR 
observed SEDs and flux levels. 

4.2. Optically thick dust clouds in the cone 

The "local" SEDs of the emission displayed in Fig.[6] pro- 
vide another key-element for understanding the MIR emission: 
the SEDs in the North-East side have a steep decline towards 
shorter wavelengths and show silicate 9.7 /im feature in ab- 
sorption. Both the SED steepness and the depth of the silicate 
feature increase with distance from the CE. On the contrary, in 



the South- West side, where the SEDs appear to be flatter, the 
9.7 yum feature is seen in emission. In other words, observations 
suggest MIR optically thick dust emission in the North-East 
cone and optically thin emission in the South- West cone. The 
most natural way to explain this behavior is through the pres- 
ence of dust within the conical regions, organized in individual 
optically thick (in the MIR) clouds, and having a small cover- 
ing factor. Indeed, in such a configuration each cloud would see 
directly the CE, and would have a much higher dust tempera- 
ture on its side facing the nucleus than on its opposite side. As 
a consequence, the overall emission from each individual cloud 
would be highly anisotropic. In the North-East cone (approach- 
ing us) we observe the clouds mostly from their colder surface, 
while the contrary holds true in the South- West cone (reced- 
ing from us). Therefore, we expect to see exactly the trend de- 
picted in Fig.0 the North-East side shows SEDs similar to that 
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Fig. 7. This figure i s a partial reproduction of Fig. 2 of 
iNenkova et al.l ( 120021) . reported here for convenience. It shows 
the predicted SED for an externally illuminated thick (jy = 
100) cloud. The varying parameter a is the angle between the 
observer and the cloud, as seen from the CE. For a = 0°, the 
cloud is located between the CE and the observer, and hence, 
the observer faces the side of the cloud which is not directly il- 
luminated by the CE. In that case, the SED is steeply decreasing 
towards shorter wavelengths and the silicate feature appears in 
absorption. For a = 180°, the CE is located between the cloud 
and the observer, and the observer faces the illuminated side 
of the cloud. The SED is much flatter than for a — 0° and the 
silicate feature appears in emission. 



labeled a = 0°, while the South- West SEDs are similar to those 
of the kind labeled a = 180°. 

A picture where the dust is organized in MIR optically thick 
clouds gives a natural explanation to both the observed emis- 
sion levels, as well as to the SED shapes, including their slope 
and silicate feature. The surface brightness of the extended 
MIR emission can be explained by this model, without invok- 
ing additional heating sources. Yet, the observations do not dis- 
card either the presence of other sources of heating, which we 
discuss in the coming section. 

4.3. Possible additional heating mechanisms 

In the former section, we discarded a picture where diffuse dust 
inside the ionizing cone, heated by the CE, was responsible for 
the observed MIR emission. On the contrary, the local SEDs 
between 7 and 20 /mi strongly favor a model in which optically 
thick dusty clouds, with small covering factor, produce the MIR 
emission. In this context, the radiation from the CE is the major 
source of heating for these clouds, but the possibility and evi- 
dence of additional heating mechanisms ought to be discussed 
here. 

First, the radiation from the CE might be beamed. This 
would add a continuum source in the direction of the beam 
enhancing the dust emission along this region. 



A second possibility would be additional heating of the dust 
by a starburst. But in the case of NGC 1068, the observed star- 
burst ring is too far away (10"-20") from the MIR structures 
discussed here. An alternative possibility could be the presence 
of local star formation along the jet (jet induced star forma- 
tion). Our superposition of the small scale MIR and radio maps 
(Fig E} shows however that the radio emission is much more 
spatially focused than the MIR emission, and consequently, the 
hypothesis of jet-induced star formation is not well supported 
by observations. 

Finally, interactions between radio and optical structures 
generate shocks that can heat the gas to very high temperatures 
(> 10 6 K). This gas then becomes an efficient source of UV and 
X-ray photons able to ionize the surrounding material. Among 
others, an expected effect is extra heating of the dust. There are 
indeed some arguments in favor of a contribution of heating by 
shocks. 

— the subarcsec NIR continuum emission jRouan et all 

Il998t iMarco & Alloini l2000h is elongated along an axis rem- 
iniscent of the radio emission axis, and particularly, the 
3.5 yum emission shows a change of direction which follows 
that of the radio jet (see Fig.@]e). 

- the sub arcsec resolu t ion Ch andra (0.25-7.50 keV) image 
published bv lYoung et alJ J200ll) reveals a bright nucleus and 
bright extended emission up to 5" to the North-East, probably 
coincident with the extended North-East radio lobe observed 
at 6cm. Moreover, this image shows that the X-ray emission 
immediately to the South of the nucleus has a more or less 
North-South extension, just like the radio emission. Could this 
X-ra y emission be contr ibuted by the shocked gas? 

- lAxon et ail Jl998h find signs of interaction between the 
radio and optical structures. For example, they find that emis- 
sion lines are split into two velocity components separated by 
~1500 kms 1 at the position of the jet axis. Outside this re- 
gion, the lines are quiescent. Supporting the shock model is 
the finding by these same authors of an excess of local optical 
continuum associated with the radio jet. The existence of an ex- 
cess of UV continuum is supported by the increase in excitation 
along the jet axis (outside the jet, the gas excitation is much 
lower). This continuum exc ess should naturally contribute to 
the dust heating, as shown bv lViHar-Martrn et alJ d200lh . Clues 
that shocks are present also come from the H? 2.12 urn an d 
Bry high resolution spectroscopy of lGalliano & Alloini J2002|) . 
At the position of the North-Bry-Fb knot, corresponding to 
the position of NLR-cloud E, the H2 emission is double 
peaked (~ 400 km s velocity difference), showing a kinemat- 
ical perturbation at this location, possibly induced by shocks. 
Bry is also very broad at this location, displaying a FWHM 
of ~1000 kms" 1 , against ~200 kms _I forH2. Can we explain 
this difference with shocks? In the absence of entrainment pro- 
cesses, the cores of the clouds impacted by the jet will be ac- 
celerated to a velocity V given by: 



where V s is the velocity propagation of the shock in the hot 
phase (that can be identified to the jet advance speed) , and «/, 
and n K are the densities of the clouds (ionized or molecular) 
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and of the hot phase (~ 10 7 K). This equation comes from ap- 
plying the balance between the pressure in the shocked inter- 
cloud medium (P/, ~ «/, x Vf), and the pressure b ehind the 
shock in the cloud (P w ~ n w x V 2 . lKlein et al 1 11994 . 

Molecular clouds will be accelerated to lower velocities 
than ionized clouds, due to the higher density (the shock finds 
higher resistance). Comparing the two velocities: 

Spectroscopy of the ionized regions (i.e. Bry emitting gas, but 
not shocked) not intera cting with the radi o jet implies densi- 
ties of a few 100 cm " 3 jAxon etall fl998h . With the molecu- 
lar gas being ~ 100 times denser, then V B r r /V H2 ~ V100 = 
10. Therefore, the ionized clouds can easily be accelerated 
to speeds ~10 times higher than the molecular clouds. This 
would explain the observed difference in the spectral profiles 
of H2 and Bry in the North-Bry-H2 knot. An alternative pos- 
sibility is that only molecular clouds impacted by slow shocks 
can survive, while those enduring high velocity shocks are de- 
stroyed. In this case, the H2 emission would come from more 
quiescent clouds. 

Therefore, it is likely that shocks, witnessed by the kine- 
matics of the molecular material, play a role in the heating of 
the dust, although this remains to be quantified. 

4.4. Conclusion 

The contribution of shocks to the heating of dust and molecu- 
lar material looks plausible and promising to explain detailed 
kinematical features. However, we find that the main source of 
extended MIR emission to the North-East and South- West are 
MIR optically thick dust clouds distributed within the ionizing 
cone, with small covering factor. 

5. The molecular emission perpendicular to the 
cone axis 

The material in the plane perpendicular to the ionization 
cone axis is visible essentially through molecular line emis- 
sion of the warm/cold material. No conspicuous NIR or MIR 
continuum emission is detected along the East-West direc- 
tion. In particular, the strong CO(2-l) and H2 2.12yum emis- 
sion lines witness the presence of m olecular material up to 
100 pc along the East- We st direction ( Schinne rer et all 12000: 
Gallia no & Altoirl [2002). The structures revealed through 
these lines correlate remarkably well. 

The 0.5" reso l ution H? 2.12yum emission line map of 
iGalliano & Alloinl J2002I) shows that the H2 emission is dis- 
tributed in two unequally bright knots, located at 1" (70 pc) 
from the CE, the brightest to the East and the weakest to the 
West. Each of the two knots displays a North-South exten- 
sion which even extends further like an arc or spiral. The East- 
H2 knot is 3 times more intense than the West-H2 one. The CO 
emission line map has a similar aspect but presents a smaller 
East-knot to West-knot intensity ratio (1.5 for CO against 3 for 
H2). The fact that the CO/H2 ratios are different in the East 
knot and the West knot suggests that extinction plays a role. 



There is another reason to believe that differential extinc- 
tion is likely to occur in NGC 1068, from the observat ional 
facts reported on h ereafter. Several studies JSchinnerer et all 
20001 iBakerl EoOO ) argued, on the basis of the CO kinematics, 
about the presence of a warp in the molecular disc surrounding 
the AGN of NGC 1068. Moreover, the maps at different scales 
(see Fig. |3 show that the PA of the main axis of the observed 
North-South structures change with scale, which supports the 
idea of a warp. X-ray observations show that the material to- 
wards the central engine is Compton thick, which implies a col- 
umn density > 10 24 crrT 2 . This Compton thick material (which 
we will call the X-ray absorber hereafter) is likely to be lo- 
cated at the very central parts of the AGN, and might be e.g. a 
corona of hot gas such as th e atmosphere of the accretion disk 
JCollin-Souffrin et all fl9961. Since the excitation of molecular 
emiss ion is dominated by the X-rays from the CE JMalonevt 
1997), a tilt between the orientations of the X-ray absorber and 
the molecular disk, which would naturally occur in the pres- 
ence of a warp, would break the axisymmetry and hence could 
significantly affect the observed molecular line intensity distri- 
bution. 

5.1. Differential extinction between the Eastern and 
Western knots 

If one assumes that both the East and West knots have in- 
trinsically the same intensity, and that the western emission 
is dimmed by dust extinction only, a column density towards 
the West-H 2 knot of the order of N Wes t = 0.5 10 22 cirT 2 (for 
A K IA V = 0.11 and A V /N H = 0.19 10~ 22 mag (nuclei cm 4 )-') 
is required. This is of the order of the column density through 
typical molecular clouds in our Galaxy. Dust extinction affects 
the H2 line photons, but do not affect the CO line photons at 
1.3 mm. Still, both CO and H2 photons can be absorbed by 
molecules in intervening clouds. Large and dense molecular 
clouds can in fact screen completely the H2 and CO photons. 
The estimation of the optical depth through an Nh = 10 22 cm -2 
cloud with «co/ w H2 = ^ 4 resu l ts m T > 10 for photons with 
v{—2crturb) < v < v(2o~tu r b), where v is the frequency of the 
crossing photon, cr tur b is the velocity dispersion of molecules 
in the cloud, and v(2cr, ur b) is the frequency of the molecular 
transition for a molecule at radial velocity of 2o- tur t- 

In the following, we show quantitatively that a small tilt 
between the molecular disc and the X-ray absorber can lead to 
a notable intensity difference in the H2 and even, but to a much 
lesser extent, in the CO emission, between the East and West 
knots. 

5.2. The model of X-ray excitation 

The excitation of the molecular emission is computed using 
the results of the model of X-ray irradiated molecular gas de- 
veloped by iMalonev et all Jl996h . This model computes the 
physical and chemical state of dense neutral gas exposed to 
an intense X-ray flux . We use the numerical predictions of 
IMalonev et all i 19961) for a power-law X-ray source of in- 
dex a = -0.7 (F Y oc y") and consider molecular clouds of 
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density n = 10 cm through which the column density is 
N = 10 22 crrT 2 . The most important parameter that controls the 
physical conditions of the molecular material is the ratio of the 
attenuated X-ray flux to density. This parameter is expressed in 
terms of an effective ionization parameter: 



&// ~ ioo 



L44 



0.9 



n 5 r l pc N 22 



where the 1-100 keV X-ray luminosity of the central source is 
Lx = 10 44 L44 erg s _1 , the distance to the X-ray source is r pc pc, 
the gas density is n — ns 10 5 cm -3 and the attenuating column 
density to the X-ray source is N n tt = \0 22 N22 cm~ 2 . The predic- 
tions of Malo nev et al.l \ 19961) for the emergent surface bright- 
ness of a cloud of density 10 5 crrT 3 in the 2.12 fim line and CO 
rotational lines (all lines together) are reproduced on Fig.|8] 



In the case of NGC 1068, L 4A = 1 (iMalonevl Il997t 
lMattetallll997l) . Since the maxima of H2 emission occur at 
l"(70pc) from the CE, we know that the value = 0.026 
(corresponding to the value of g e ff at the maximum of I/N22 
for H2 in Fig.|8j corresponds to r pc —10. This gives us access to 
the column density between the CE and the H2 knots, using the 
expression of ^ e ff- 



N22 



IOOL44 

Z eff n 5 d 2 ) 



0.75 



Therefore, this value is much smaller than the column density 
of > 10 24 cirr 2 to the CE, deduced from X-ray observations. 
This supports the hypothesis that the central X-ray absorber 
and the emitting molecular material are not coplanar. Hence, 
we have developed a new model along this line of arguments, 
which is presented in the next section. 




10° 10 1 10 2 10 3 10 4 10 



Fig. 8. Reproduction of the results for the emergent intensity in 
the H2 2.12 fim emission line and t he CO rotational lines for a 
cloud of column density 10 22 cirT 2 jMalonev et allll996l) . See 
Sect.E2]for details. 



to the observer 

non-obscured 
emitting region 




central Xray 
absorber 

Fig. 9. Sketch of a configuration consisting of a molecular disc 
(grey area) and an X-ray absorber (black cone) tilted with re- 
spect to the disc. See text for details. 



5.3. Molecular disc model 

The configuration sketched in Fig. [5] consists of a molecular 
disc and of a central axisymmetrical X-ray absorber, tilted with 
respect to the molecular disc. In Fig. [9] the grey rectangle rep- 
resents a cut through the edge-on view of the molecular disc. 
The central X-ray absorber axis is tilted with respect to the axis 
of the disc. The region filled with the "hexagon" patterns are 
the regions where the X-rays are deeply attenuated, and con- 
sequently from which no strong emission can emerge. The cir- 
cle centered on the CE traces the radius at which, in the ab- 
sence of a central X-ray absorber, the value of g e ff would cor- 
respond to the maximum of emission. The regions from which 
the emission will actually arise are shown with the hatched cir- 
cles. The observer sees directly the "non obscured emitting re- 
gion", while it sees the "obscured emitting region" through the 
molecular disc. This sketch provides a qualitative explanation 
to the fact that the molecular emission is offset from the CE, 
and how it can naturally appear brighter on one side. 

We have built a simple numerical model to simulate such a 
configuration. It consists of the following two components: 

(1) A disc distribution of molecular clouds around the cen- 
tral X-ray source. All clouds are identical to the densest cloud 
studied bv lMalonev et alJ {l996), i.e. their density is 10 5 cm -3 
and their column density is N = 10 22 crrT 2 . Complex cloud 
number density distributions describing flared or warped discs 
were computed as well, but a simple distribution (with less free 
parameters) is sufficient for getting interesting results. In an 
(x,y,z) reference frame, where z is the axis of the disc, we have 
chosen the following distribution: 

d c (z)=D c e-^ A 

where d c is the cloud number density, D c is the density at z=0, 
and (Z, A) controls the dependence of the density on the height 
above the equator of the disc. A flat rotation curve was as- 
signed to the disc (plateau at 130 kms -1 ), in order to estimate 
the opacity of the clouds in the lines. 

(2) A central axisymmetrical X-ray absorber, for which the 
equatorial attenuation is very high (> 10 24 cirT 2 ). We have set 



E. Galliano et al.: Location and environment of the central engine in NGC 1068 



13 



clump density [pc 3 ] l°g(£err) 




(c) -2 -1 1 2 (d) -2 -1 1 2 



Fig. 10. (a) cloud number density distribution for an edge-on view. The image represents the density in the plane perpendicular 
to the disc and passing through the center. The unit is cloud per pc 3 . (b) The values of g e ff for the same plane as in (a). The 
X-ray absorber is slightly rotated with respect to the disc, (c) Final orientation of the molecular disc, (d) Final orientation of the 
X-ray absorber. 



the azimuthal dependence of its attenuating column density to 
the following form: 

AW©) = AUO) 10 rf(Q) , 

where © is the azimuthal angle and N alr (Q) the equatorial thick- 
ness of the X-ray absorber, that is N att > 10 24 cm -2 . 

These two components can be freely and independently ro- 
tated in space. 

For each cell of a 3D grid we have computed the value of 
the attenuation towards the X-ray source, as the sum of the at- 
tenuations due to both components. The value of g e ff in each 
cell is then computed. 

The attenuation of the emission line photons on their way 
to the observer is computed taking into account two effects. (1) 
the dust extinction for H2 only, (2) the high opacity of some 
intervening molecular clouds in both H2 and CO lines (see 
AppendixIXl. 

5.4. Results 

Fig.[IT5]and Fig.^Jdisplay the results of an illustrating model 
for a seeing of 0.4". 



In Fig. [TO] the top two figures (a) and (b) present an edge-on 
view of the disc. For a plane perpendicular to the disc, passing 
through the center, image (a) shows the cloud density distri- 
bution, and image (b) the values of for an X-ray absorber 
tilted with respect to the molecular disc by 15° in the direction 
of the x axis. The subsequent simulated images (c,d,f) are not 
seen edge-on anymore. Images (c) and (d) respectively show 
the actual orientations of the chosen model for the molecular 
disc and the X-ray absorber. 

Fig. [^compares model predictions with observations: (a) 
displays the H2 emission map and (b) the model prediction for 
the H2 emission, (c) displays the CO(2-l) emission map and 
(d) the model prediction for the CO (all lines) emission. 

The parameter values corresponding to this model are: 
D c = lOOpc" 3 , Z = 20pc, A = 1.15, N all (0) = 10 26 ctrT 2 , 
(3 - 1.5, cr cin - 100 kms~' and cr Hlrb - 20 kms -1 

The orientation of the two components is illustrated by the 
ellipses of Fig.^3(c) an d (d). 

This two component model reproduces quite well the most 
characteristic features of the observations of the following 
H2 2.12yum and CO emission: The distance of the observed 
emission knots to the CE, the extended arc or spiral structure, 
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Fig. 11. (a) Observed H2 surface brightness. The unit is 10 3 erg s _1 cuT 2 sr . (b) Modeled H2 surface brightness, (c) Observed 
CO(2-l) surface brightness (integrated between -100 and 100 km s _1 ). (d) Modeled CO rotation lines surface brightness. The unit 
is 10" 3 erg s _1 crrT 2 sr -1 . In all figures, the x and y axis unit is arcsec. 



the intensity ratios between the East and West knots in H2 and 
CO (3.2 for H2 and 1.25 for CO) and is consistent with the 
Compton thickness of the absorber along the line-of-sight. Yet, 
the surface brightnesses of the simulated emission knots are too 
large compared to the observations. Choosing a smaller cloud 
number density lowers the intensity, but lowers as well the in- 
tensity difference between the East and West knots. This sug- 
gests that there is probably a difference in the amount of mate- 
rial between the East and the West knots. If we consider such a 
difference, we can then fit more easily the surface brightnesses 
and the intensity (all lines together) ratios at the same time. 
Still, in that case a configuration with selective extinction is 
required. 

Another point is that, for the chosen model, the predicted 
CO emission knots are more extended than observed. This can 
be resolved by adjusting the dependance of the cloud density 
distribution on the distance to the central engine. The fact that 
we predicted the total CO rotational emission intensity explains 
part of this discrepancy. We also point out that, because CO is 
still emitted at lower g e ff values than H2, in the configuration 
we present, the H2 and the CO peaks do not fall exactly at the 
same position, the CO emission being slightly rotated counter- 
clockwise with respect to the H2 emission. 



An interesting observational feature, which is not repro- 
duced by this model, is the increase of CO emission towards 
the South in the western knot, a feature which is not observed 
fo r the H? emission. In the k inematical configuration proposed 
bv ldalliano & Alloinl J2002ll to fit the H2 line profiles, both a 
rotating and an outflowing kinematical components had to be 
introduced. One notices that it is precisely in the South- West 
and North-East regions that the projected velocities of both 
components have the greatest difference (where the profiles are 
the most asymmetric). In the western region, the outflowing 
component would be more absorbed by the rotating disc to the 
North (the projected velocities tend to be more equal), than to 
the South. Hence, this would qualitatively explain a rise of CO 
emission to the South of the western knot. 

6. Summary and conclusion 

We have discussed in this paper several features of the central 
100 pc region around the CE in the AGN of NGC 1068. 

The second section was dedicated to revisiting the relative 
positioning of maps at different wavelengths. We provide a new 
estimate of the absolute positioning of the unresolved peak in 
the K-band. After compiling the most important observations 
published, and after elucidating the apparent disagreement be- 
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tween different results, we could register with a high precision 
the radio, NIR, optical and UV continuum maps, as well as the 
submillimiter, NIR, optical and UV emission line maps. We 
produced a series of map overlays, which constitute an essential 
tool for giving a general view of the complexity of the central 
region of NGC 1068. These superpositions show that the MIR 
is very well correlated with the [OIII] emitting NLR clouds, 
revealing the survival of dust grains in the UV-irradiated ion- 
ization cone. The final registration of the small scale high res- 
olution radio maps with respect to the HST optical maps show 
a surprising anti-correlation: no NLR cloud is found along the 
radio jet trajectory. This suggests a picture where the radio jet 
sweeps away material within and along the axis of the ionizing 
cone. 

The third section deals with the modeling of the NIR-MIR 
unresolved peak emission in NGC 1068. First, we could show 
that the SED fitting approach is somewhat limited since several 
acceptable combinations of torus geometries and dust parame- 
ters are compatible with the observed IR SED of NGC 1068. 
We conclude that no strong constraint on the geometry of 
the torus can be set using this approach, while its extension 
and thickness can to some extent be estimated: the torus in 
NGC 1068 is smaller than 10-20pc and moderately thick in the 
MIR. 

Then, a discussion is provided about the origin of the ob- 
served extended MIR emission, mainly along the axis of the 
ionizing cone. We find that the most natural way to match the 
observations is to consider individual optically thick clouds, 
distributed in the NLR and having a small covering factor. An 
alternative picture, consisting in diffuse dust distributed within 
the cone, is confidently discarded. Possible additional sources 
of heating are discussed, and we show in particular that the 
likely presence of shocks can explain some of the observed 
kinematical features such as the difference between the FWHM 
of the molecular emission lines and of the ionized gas emission 
lines, in the North-Bry-H2 knot in particular. 

The last section presents a model of the excitation and dis- 
tribution of the H2 and CO line emission along the East- West 
direction (perpendicular to the ionizing cone axis). We argue 
that the observed intensity differences between the eastern and 
the western regions, with respect to the CE, can be partly due to 
selective extinction produced by a slight tilt between the orien- 
tations of the Compton thick X-ray absorber in the very central 
region of the AGN and the molecular material distribution at 
larger scales. 
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Pi = Seff/Sceii- In Plxl00% of the cases, cloud/, can be con- 
sidered to be hiding completely cloud a , and in (1-Pl)xl00% of 
the cases, cloudb does not hide cloud a at all. 

We can define as well V e /f (the analogous to r e ff in the 
velocity space) such that, if -V e ff ^ V a -Vi, < V e ff, we can 
consider that the photon, in average will be absorbed, where V a 
and Vt are the respective velocities of cloud a and cloudb along 
the line of sight. We find V e /f = \ .2o-, ur h- 

Let o- C i n be the velocity dispersion between the cloud ve- 
locities. If AV = absiVa - Vb) is the difference between the 
radial velocities of cloud a and cloudb, then the probability for 
cloudb to have the same velocity as cloud a is: P2 = V e ff(o- cin - 

Finally, in (PI x P2)% of the cases, the photon emitted by 
cloud a will be absorbed by cloudb- Thus, we approximate the 
attenuation by molecules in intervening clouds, by a factor (1 - 
Pi x Pi) per cloud in the cells located along the line of sight to 
the observer. 



Appendix A: Attenuation of H 2 and CO photons by 
intervening molecular clouds 

We computed that a molecular cloud of column density N = 
10 22 cirT 2 is completely opaque (t > 10) to CO and H2 pho- 
tons, if the radial velocity difference between the emitting 
molecule and the molecular cloud is < 2o- tur b, where cr tU rb is 
the velocity dispersion inside the cloud (of the order of 10- 
20kms~'). Let us estimate the corresponding attenuation ef- 
fect. 

We consider two cells of the grid aligned along the line of 
sight, cell a and cellb- Cell a contains one emitting cloud called 
cloud a and cellb contains one non-emitting cloud called cloudb- 
Cellb is located between cell a and the observer. We can esti- 
mate a cross section S e ff of radius r e ff for cloud a , so that we 
can consider that, when the center of cloudb lies inside S e /f, 
cloud a is completely covered, and when the center of cloudb 
lies outside, then cloud a is completely uncovered. The quantity 
r e ff is such that the overestimation in covering when the center 
of cloud/, does not coincide with the center of cloud a (but still 
lies within S e /f) is balanced by the underestimation in cover- 
ing when the center of cloudb is outside S e /f (but still covers 
part of clouda). We find r e //=0.89r e / OM j. 

Then we can simply compute the absorption in terms of 
probability for the center of cloudb to lie within the cross sec- 
tion of cloud a in the area of a cell S ce u- This prbability is 
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